clear

use "1920_1940_merged.dta"

set seed 349132

preserve
keep if age1920>=6 & male1920==1
keep if (white1920==1|black1920==1)

*Drop if dad's location not found
keep if occ1950_pop1920~=.

*Drop if child's occupation is missing, or if father's occupation is missing
drop if occ19501940==999
drop if occ19501940==997
drop if occ19501940==979

drop if occ1950_pop1920==999
drop if occ1950_pop1920==997
drop if occ1950_pop1920==979

tab age1920
tab age1940
count

egen taghistid = tag(histid1920)

sort age1920 taghistid
by age1920 taghistid: gen scale = _N if taghistid==1
by age1920 taghistid: egen sonrankx = rank(occscore1940) if taghistid==1
by age1920 taghistid: egen dadrankx = rank(occscore_pop) if taghistid==1
replace sonrankx = sonrankx/scale
replace dadrankx = dadrankx/scale

bysort histid1920: egen sonrank = max(sonrankx)
bysort histid1920: egen dadrank = max(dadrankx)

egen index = group(cty_fips1920)
su index
local indexmax = r(max)
gen constantreg = 0 
gen coefficient = 0 
forvalues i = 1(1)`indexmax'{
reghdfe sonrank dadrank if index==`i' [pw=weight]
mat C = e(b)
replace constantreg = C[1,2] if index==`i'
replace coefficient = C[1,1] if index==`i'
}

gen p25 = constantreg + coefficient*.25
gen p50 = constantreg + coefficient*.5
gen p75 = constantreg + coefficient*.75

#delimit ;
local vars constant;
foreach x of local vars{;
gen `x'wt = `x'*weight;
drop `x';
rename `x'wt `x';
};

#delimit ;
collapse (mean) statefip1920 hasnormalschool p25 p50 p75 (sum) constant, by(cty_fips1920 white1920);

#delimit cr
gen constantrd = round(constant,1)
gen below10 = constantrd<10
count if below10==1
bysort statefip hasnormalschool white1920: egen totctytype = count(hasnormalschool)
bysort statefip hasnormalschool white1920: egen totbelow10 = total(below10)
gen fracbelow10 = totbelow10/totctytype
bysort statefip white1920: egen maxfracbelow10 = max(fracbelow10)

tab statefip white1920 if maxfracbelow10>=.5
drop if maxfracbelow10>=.5


#delimit cr
sort cty_fips

local vars p25 p50 p75
foreach x of local vars{

reghdfe `x' hasnormalschool if white1920==1, absorb(statefip1920) cluster(statefip1920)
su `x' if e(sample)==1 & hasnormalschool==0
local mean = r(mean)
tab statefip1920 if e(sample)==1
outreg2 using censustree19201940mearank.xls, append excel dec(3) addstat(dvmean, `mean') adec(3)

reghdfe `x' hasnormalschool if white1920==0, absorb(statefip1920) vce(robust)
su `x' if e(sample)==1 & hasnormalschool==0
local mean = r(mean)
tab statefip1920 if e(sample)==1
outreg2 using censustree19201940mearank.xls, append excel dec(3) addstat(dvmean, `mean') adec(3)

#delimit ;
permute hasnormalschool _b[hasnormalschool], strata(statefip1920) reps(1000): reg `x' hasnormalschool i.statefip1920 if white1920==0;
#delimit cr
}
